Load packages
Step 1 Load data from 10X folder
# Load in the RNA UMI matrix
cbmc.rna <- as.sparse(x = read.csv(file = "../../Tutorial/data/GSE100866_CBMC_8K_13AB_10X-RNA_umi.csv.gz",
sep = ",", header = TRUE, row.names = 1))
cbmc.rna <- CollapseSpeciesExpressionMatrix(object = cbmc.rna)
# Load in the ADT UMI matrix
cbmc.adt <- as.sparse(x = read.csv(file = "../../Tutorial/data/GSE100866_CBMC_8K_13AB_10X-ADT_umi.csv.gz",
sep = ",", header = TRUE, row.names = 1))
cbmc.adt <- cbmc.adt[setdiff(x = rownames(x = cbmc.adt), y = c("CCR5", "CCR7",
"CD10")), ]Step 2 Create object
Step 3 Pre-process
1) Filter out unwanted cells (optional)
for this dataset, we don’t need to filter out unwanted cells
2) Remove unwanted genes (optional)
for this dataset, we don’t need to filter out unwanted genes
3) Normalization
data Normalization for both ADT (CLR) and RNA (log)
4) Indentify HVGs for RNA data
Call seurat function to identify highly variable genes (HVG) for RNA data
Step 4 Linear dimension reduction (PCA)
directly call Seurat function for linear dimension reduction (PCA)
Step 5 Determine number of PCs
call Seurat function JackStraw to determine number of PCs
# cbmc <- JackStraw(cbmc, num.replicate = 100) cbmc <- ScoreJackStraw(cbmc,
# dims = 1:20) JackStrawPlot(cbmc, dims = 1:15)
ElbowPlot(cbmc)Step 6 Distance calculation and joint distance calculation
calculate cell-cell distances for RNA, ADT and joint. alpha was set to 0.5 as initial, number of PC was set to 20 by default.
Step 7 Non-linear dimension reduction (UMAP and t-SNE)
run UMAP as Non-linear dimension reduction for RNA, ADT and joint analysis.
Step 8 Clustering
Step 9 Visualization ADT vs RNA vs Joint
1) Cell clusters
# gridDimPlot(cbmc, wide.rel = 1.5, legend = FALSE, reduction.prefix =
# 'tsne_', height.rel = 0.5)
plots <- generateGridDimPlot(cbmc, legend = FALSE, reduction.prefix = "tsne_",
darkTheme = FALSE)
listPlot(object = plots)###### user also can only plot some of those plots by index, figure ident or
###### figure map info listPlot(object = plots, fig.ident = 'RNA')
###### listPlot(object = plots, fig.ident = 'RNA', fig.map = 'RNA')
###### user can use plotInfo() function to get index, figure ident and figure map
###### information, then plot figures by index
plotInfo(plots)
# listPlot(object = plots, fig.id = 1)2) Heat maps
# Heatmap for joint clusters
heatMapPlot(object = cbmc, group.by = "jointClusterID", height.rel = 1, adt.label = TRUE)# Heatmap for RNA clusters
heatMapPlot(object = cbmc, group.by = "rnaClusterID", height.rel = 1, adt.label = TRUE)# Heatmap for ADT clusters
heatMapPlot(object = cbmc, group.by = "adtClusterID", height.rel = 1, adt.label = TRUE)Step 10 Build trajectory using MST method
Step 10 MST Trajectory Visualization
p1 <- trajectoryMSTPlot(object = cbmc, group.by = "jointClusterID", reduction = "tsne_joint",
pg = "jointPG")
p2 <- pseudoTimeMSTPlot(object = cbmc, reduction = "tsne_joint", pg = "jointPG",
pseudotime = "jointTreeMSTTime")
plot_grid(p1, p2, ncol = 2)Step 11 build trajectory using kNN method
Step 11kNN Trajectory Visualization
p1 <- trajectoryKNNPlot(object = cbmc, group.by = "jointClusterID", reduction = "tsne_joint",
graph = "jointGraph", line.size = 0.1)
p2 <- clusterConnectionPlot(object = cbmc, group.by = "jointClusterID", reduction = "tsne_joint",
graph = "jointGraph", line.color = "black", cutoff = 0)
p3 <- pseudoTimeKNNPlot(object = cbmc, reduction = "tsne_joint", graph = "jointGraph",
pseudotime = "jointGraphKNNTime", line.size = 0.1)
plot_grid(p1, p2, p3, ncol = 2)